# TUTORIAL 6 : Free Energy (FE) and replica-exchange methods¶

Author: | Andrey V. Brukhno (andrey.brukhno{at}stfc.ac.uk) |
---|

## General concept¶

Most often it is impossible to calculate the absolute free energy value (we will see why later on).
However, *relative* free energy can always be expressed via the probability of the associated state:

then for the **FE difference** we get

This is the most useful statistical-mechanical approach to FE evaluation, and all *advanced* FE methods are based on the above formula (althogh there exist also *thermodynamic* approaches, usually more demanding and time-consuming).

## Example: Potential of mean force as FED¶

Potential of mean force (PMF or POMF) is defined as the work done upon reversibly bringing two particles from the state of zero-interaction into the state of full-interaction at a given distance, *W(R)*.

PMF can be directly obtained from RDF:

where the zero-interaction state corresponds to infinite separation \((R\to\infty)\).

Alternatively, one can calculte

- the work to bring two particles \((i,j)\) from \(|\Delta\mathbf r_{ij}|=R_1\) to \(R_2\).

Clearly, if \(g(R_1)=1\), we get back the total work \(W(R=\infty\to R_2)\).

### Why PMF is a free energy?¶

First, note that

We can rewrite the canonical configuration integral (configurational part of the partition function)

as an integral (i.e. sum) over the substates corresponding to \(|\Delta\mathbf r_{ij}|=R\).
This implies that we can consider each substate as an individual thermodynamic state, whereas \(R\)
serves as an additional **external constraint**, similar to *N,V,T* or *p*. Thereby,

where

Therefore, we can evaluate the FED between any two substates (different particle separations) by performing a Metropolis
sampling along *R* and collecting a histogram of hits (or visits) on a discrete grid with sufficiently small bin size,
\(\Delta R\), which is of course equivalent to calculating \(g(R)\)!

NOTE: One does not need to know the total partition function (nor the total configuration integral) as it only serves as a normalization constant which cancels out!

- Question: what is the
*“unknown constant”*in the first equation above? - Similar formulae hold true generally, for other constrained substates (see the Appendix below).

## FED methods in DL_MONTE¶

In general, Metropolis sampling is prone to getting trapped (stuck) in local FE minima and, hence, suffering from undersampling (even avoiding!) regions around energetic or entropic barriers. This makes estimating FED via probablity ratio troublesome…

The solution to this problem is to **bias the sampling** away from the minima and towards the barriers.

**Ideally the bias potential is to perfectly compensate for the underlying free energy differences!****This implies that the resulting probablity distribution in a perfectly biased simulation is uniform, i.e. flat!**

## Parallel tempering with replica-exchange (RE)¶

**Two approaches to parallelization in DL_MONTE**

In a **parallel multicanonical simualtion** each workgroup treats the system at a different tempreature
and generates it own MC trajectory. It is then possible to arrange periodic exchnage of configurations between
the workgroups, which is equivalent to concurrent jumps in *T* with the acceptance probability:

If a FED evaluation is carried out in the same simulation, say FED(R), the acceptance rule includes an extra term:

## Exercises¶

Now try the exercises to learn how to perform FED calculations with DL_MONTE:

TUTORIAL 6 : Exercise 1 - Reproducing interaction potentials

U(r)with EE & WL/RE methodsTUTORIAL 6 : Exercise 2 - Calculating PMF for two charged colloids in ionic cloud

TUTORIAL 6 : Exercise 3 - Umbrella sampling in subranges (windows) with the use of WHAM method

## Appendix: FED via a **Generalised (Expanded) Ensemble**¶

The biasing is best illustrated with a *generalised, so-called Expanded Ensemble (EE) approach*.

The EE is defined by a generalised partition sum:

or in a discretised form:

where \(\{{\mathbf q}\}\) stand for one or more **generalised coordinates (or constraints)**, which can also include normal thermodynamic variables (e.g. *N,V,T, p* etc), and \(w({\mathbf q})\) is a **biasing weight** assigned to a particular
set of \({\mathbf q}\).

Clearly, if Metropolis sampling is carried out in such a generalised (so-called expanded) ensemble, including MC moves in \({\mathbf q}\), the probablity of an EE substate

and the *effective* (biased) free energy difference

Finally, if the weights \(\{w({\mathbf q}_m)\}\) are chosen so that \(p({\mathbf q}_m)/p({\mathbf q}_n) = 1\), we get:

**Ideally the EE weights are to perfectly compensate for the unbiased probablity distribution!**

In practice, though, one has to arrange an iteration to self-consistently update the weights until all the substates are populated reasonably well, after which the FED is obtained as

NOTE: It is convient (and conventional) to define the EE weights via a *biasing potential*:

so that

**Ideally the bias potential is to perfectly compensate for the underlying free energy differences!**